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Annual  Summary  Report 

Quantum  Monte  Carlo  for  Molecules 

Principal  Investigators: 

William  A.  Lester,  Jr. 

Peter  J.  Reynolds 

Materials  and  Molecular  Research  Division 
Lawrence  Berkeley  Laboratory 
University  of  California 
Berkeley,  California  94720 

Description  of  Problem  and  Approach 

Introduction. 

In  the  last  few  years,  Monte  Carlo  approaches  have  enjoyed  ever  wider 
application  in  the  treatment  of  problems  in  theoretical  physics  [l].  The  Monte 
Carlo  method  is  statistical  in  nature,  based  on  the  generation  of  “random” 
numbers  or  “coin  tosses.”  Thus  it  is  perhaps  easy  to  imagine  using  the  Monte 
Carlo  method  for  treating  inherently  statistical  models.  It  is,  however,  also 
readily  applied  to  problems  of  numerical  integration  [2].  How  to  solve  m any¬ 
body  equations  such  as  the  Schrodinger  equation  by  Monte  Carlo  is  far  less  obvi¬ 
ous.  Nevertheless,  many-body  problems  and  the  Schrodinger  equation  are  readily 
treated  by  Monte  Carlo. 

Our  focus  here  is  on  the  quantum  mechanical  Monte  Carlo  (or  OMC) 
methods  (3).  These  are  the  methods  that  explicitly  treat  the  quantum  aspects  of 
a  problem,  such  ;is  solving  the  Schrodinger  equation.  'Phis  is  done  statistically  by 
the  simulation  of  an  appropriate  random  process.  In  particular,  the  formal  simi¬ 
larity  of  the  Schrodinger  equation  with  a  ”  fusion  equation  allows  one  to  calcu¬ 
late  quantum  mechanical  expectation  values  as  Monte  Carlo  averages  over  an 
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ensemble  of  appropriately  chosen  random  walks.  The  method  is  proceduraily 
quite  simple.  As  a  result,  QMC  provides  an  attractive  alternative  to  the  conven¬ 
tional  variational  and  perturbation-theoretic  techniques  used  in  physics  and 
chemistry. 

Past  Work 

We  have  applied  QMC  sucessfully  to  the  calculation  of  the  total  energy  of  a 
number  of  molecular  systems  [4-7].  In  every  case  we  have  achieved  very  high 
accuracy  compared  with  experimentally  inferred  and  exact  results  (where  avail¬ 
able),  as  well  as  ab  initio  configuration  interaction  calculations.  In  most  cases, 
OO-lOO/o  of  the  correlation  energy  has  been  obtained. 

More  recently  this  research  program  has  also  begun  to  make  a  significant 
contribution  [5-7]  to  obtaining  quantities  of  such  physical  interest  as  electron 
affinities,  ionization  potentials,  binding  energies,  reaction  barriers,  level  splittings, 
and  so  on.  Since  much  of  chemistry  takes  place  predominantly  in  the  valence 
electrons  of  a  system,  the  quantities  of  interest  are  usually  small  differences  of 
large  total  energies.  Calculations  for  such  difference  quantities  are  a  far  more 
difficult  task  for  Monte  Carlo,  since  a  small  statistical  uncertainty  (e.g.  of  as  little 
as  0.1. o)  in  the  separate  total  energies  can  mask  the  sought-after  energy 
difference.  To  reduce  the  statistical  error  to  the  level  needed  by  "brute  force"  is 
costly  in  computer  time,  as  the  standard  deviation  decreases  only  as 
{('Pf  tune)*.  Algorithmic  developments,  such  as  differential  QMC  [S]  have 
undergone  preliminary  development  and  hold  promise  for  reductions  in  variance 
through  correlated  sampling  techniques.  Another  avenue  utilizes  the  zero- 
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variance  property  of  the  approach:  the  variance  vanishes  as  a  trial  wave  function 
'I' T  approaches  the  true  eigenfunction.  To  attempt  to  take  advantage  of  this,  we 
have  previously  developed  an  iterative  procedure  for  improving  ^  T  [7]. 


t 

j 

I 

I 

t 

I 
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Recently,  we  have  also  devised  improved  forms  of  correlation  factors  for  use  with 
SCF  trial  wave  functions,  and  methods  of  Monte  Carlo  optimization  of  the  corre¬ 
lation  parameters. 

More  specific  examples  of  systems  we  have  treated  include  selected  points 
along  the  reaction  coordinate  of  the  H  +  H2  exchange  reaction  [6],  with  particular 
emphasis  on  the  saddle-point  geometry  for  which  Liu  [9]  has  performed  the  most 
extensive  Cl  calculation  to  date.  The  bound  for  the  barrier  height  which  we 
obtained  by  QMC  is  0.16  kcal/mole  below  Liu’s  bound,  and  probably  lies  within 
0.1  kcal/mole  of  the  exact  answer.  In  addition,  we  were  able  to  obtain  these 
results  with  only  single-determinant  trial  functions,  and  a  basis  set  expansion  at 
only  the  double-zeta  level. 

I’ntil  recently,  most  QMC  applications  have  been  limited  to  calculate  of 
energies  of  ground  states  and  lowest  states  of  a  given  symmetry.  As  an  example 
of  the  latter,  we  calculated  the  singlet-triplet  splitting  in  methylene  [7],  The  two 
states  studied  are  both  lowest-energy  states  of  their  respective  symmetries.  In 
these  studies,  the  accuracies  obtained  with  QMC  have  been  comparable  to  the 


hest  achieved  by  conventional  ab  initio  methods.  In  the  next  section,  we 
our  recent  work  extending  QMC  to  excited  states  of  the  same  symmetry 
ground  state.  We  also  report  on  calculations  of  properties  unrelated 
energy  such  as  moments  of  the  charge  distribution,  and  transition 


discuss 
as  t  he 
to  the 
dipole 
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moments. 

In  addition  to  properties,  it  is  certainly  useful  to  be  able  to  compute  the  ana¬ 
lytic  gradient  of  the  energy  with  respect  to  nuclear  coordinates  in  order  to  deter¬ 
mine  forces,  and  thereby  equilibrium  and  transition  state  geometries  [10]  and  (by 
finite  difference  or  higher  analytic  derivatives)  harmonic  vibrational  frequencies 
[l  1  ] .  In  the  recent  past  we  have  implemented  QMC  methods  for  doing  these  cal¬ 
culations.  They  are  currently  working,  but  need  improvement,  especially  with 
respect  to  efficiency.  Using  this  approach,  we  have  performed  calculations  on  H9 
at  a  few  nuclear  separations  [12].  QMC  is  competitive  with  Cl,  and  far  superior 
to  Hartree-Fock.  The  combination  of  QMC  energy  points  with  the  additional 
information  available  from  the  derivatives  has  allowed  us  to  reconstruct  the 
potential-energy  curve  of  H0  from  just  four  QMC  data  points.  The  Monte  Carlo 
potential-energy  curve  so  obtained  is  indistinguishable  from  the  exact  one  [13]. 

QMC  Theory 

Briefly,  one  simulates  a  quantum  molecular  system  by  allowing  it  to  evolve 
under  the  time-dependent  Schrodinger  equation  in  imaginary  time.  It  is  easy  to 
show  [  J]  that  the  use  of  imaginary  time  causes  the  system  to  approach  a  station¬ 
ary  state  which  is  the  lowest  state  of  a  given  symmetry.  Many  properties  may 


then  be  ‘‘measured'’  as  averages 


over  the  resulting  equilibrium  distribution. 


By  writing  the  imaginary-time  Schrodinger  equation  with  a  shift  in  the  zerc 


of  energy  ?us 


r  J  v|/  ( [{  (  ) 

- =  D  ./  )  +  \Et  l '(£)]*(£,/ 


/y. 


we  see  that  this  equation  is  a  generalized  diffusion  equation.  The  first  term  on 
the  right-hand-side  is  the  ordinary  diffusion  term,  while  the  second  term  is  a 
position-dependent  rate  (or  branching)  term.  For  an  electronic  system, 
D  =/i2/2mf  ,  R  is  the  three-N  dimensional  coordinate  vector  of  the  N  electrons, 
and  V  ( R  )  is  the  Coulomb  potential.  Since  diffusion  is  the  continuum  limit  of  a 
random  walk,  one  may  simulate  Eq.  (l)  with  the  function  'k  (note,  not  vk2)  as  the 
density  of  “walks”.  The  walks  undergo  an  exponential  birth  and  death  as  given 
by  the  rate  term.  This  connection  between  a  quantum  system  and  a  random 
walk  was  first  noted  by  Metropolis,  who  attributes  it  to  Fermi  [14]. 


The  steady-state  solution  to  Eq.  (1)  is  the  time-independent  Schrodinger 
equation.  Thus  we  have  '&{R  ,t  )—*<t>(R),  where  <$>  is  an  energy  eigenstate.  The 
value  of  Et  at  which  the  population  of  walkers  is  asymptotically  constant  gives 
the  energy  eigenvalue.  Early  calculations  employing  Eq.  (I)  in  this  way  were 
done  by  Anderson  on  a  number  of  one-  to  four-electron  systems  [15]. 

Unfortunately,  in  order  to  treat  systems  larger  than  two  electrons,  the  Fermi 
nature  of  the  electrons  must  be  taken  into  account.  The  antisymmetry  of  the 
eigenfunction  implies  that  *  must  change  sign;  however,  a  density  (e.g.  of  walk¬ 
ers)  cannot  be  negative.  To  handle  this,  Anderson  made  simplifying  assumptions 
about  the  positions  of  the  nodes.  His  method  was  ad  hoc  ,  and  not  readily  gen- 
eralizable.  Another  method  which  imposes  the  antisymmetry,  and  at  the  same 
time  provides  more  efficient  sampling  (thereby  reducing  the  statistical  "noise"),  is 
importance  sampling  with  an  antisymmetric  trial  function  ^  T  (see  e.g.  Ref.  4). 

I  he  zeroes  (nodes)  of  ^  T  become  absorbing  boundaries  for  the  diffusion  process; 
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this  maintains  the  antisymmetry.  The  additional  boundary  condition  that  ^ 
vanish  at  the  nodes  of  V  T  is  the  fixed-node  approximation  [4,16].  The  magni¬ 
tude  of  the  error  thus  introduced  depends  on  the  quality  of  the  nodes  of  'I'  T{R ), 
and  vanishes  as  ^  T  approaches  the  true  eigenfunction.  To  the  extent  that  T  is 
a  good  approximation  of  the  wave  function,  the  true  eigenfunction  is  almost  cer¬ 
tainly  quite  small  near  the  nodes  of  ^  T  .  Thus  one  expects  the  fixed-node  error 
to  be  small  for  reasonable  choices  of  <if  T  .  Work  on  a  number  of  systems  has 
borne  this  out  [4-7,17,18].  In  addition,  this  error  is  variationally  bounded. 

To  implement  importance  sampling  and  the  fixed-node  approximation,  Eq. 
(1)  is  multiplied  by  'l' T ,  and  rewritten  in  terms  of  the  new  probability  density 
/  (K  ,0=  ^  r  (B-  )V{&  )•  The  resultant  equation  for  /  (R  ,t  )  may  be  written 

-§7  -  D  +{ET-EL{£)\f  -DV-lfF(Z)}  .  (2) 

The  local  energy  Ei  (R_ ),  and  the  local  velocity  field  F  (R. )  are  simple  functions 
of  'k  T  given  by 

El(R)=H*t(R  )/*t(R)  ,  (3a) 

and 

F(R  )=2V*t{R)/*t(R  )  .  (31.) 

Equation  (2),  like  Eq.  (1),  is  a  generalized  diffusion  equation,  though  now  with  the 

addition  of  a  drift  term  due  to  the  presence  of  F  . 

In  order  to  perform  the  random  walk  implied  by  Eq.  (2)  we  use  a  short-time 
approximation  to  the  Cireen  s  function.  The  Cireen’s  function  is  used  to  evolve 
the  distribution  forward  in  time.  i.e.  f  (R  .t  )—/(/£'./  +  r).  This  process  is 
iterated  to  large  t  .  Such  an  approach  becomes  exact  in  the  limit  of  vanishing 


A 


*  t 

> 

*■ 


time-step  size,  r.  Asymptotically,  f{R,t )— * ■ !  oo{R  )—^  t  {B.  )• 

The  function  4>{R)  is  the  lowest-energy  eigenfunction  of  the  Schrodinger 
equation  for  the  imposed  set  of  nodes.  Although  neither  this  function  nor  } ^  is 
known  analytically,  we  can  nevertheless  sample  desired  quantities  from  the  equili¬ 
brium  distribution.  Averages  taken  with  respect  to  the  distribution  f ^  are 
known  as  mixed  averages.  For  example,  sampling  a  quantity  A  in  equilibrium 
after  N  samples  gives  the  average  (in  the  limit  of  large  X ) 


configs 


=  <A  > 


—  ffooiS. )  A  dR 


_ _/ *  T{R)d>{R)A  dR 

/*r(£  )kR)dR  ' 

or  in  abbreviated  Dirac  notation  (with  the  normalization  absorbed), 

<A  > /«,*“<  \A\i».  (5) 

On  the  other  hand,  the  correct  expectation  value  of  A  ,  for  a  state  d>,  is 

<d>  |  A  |  <?>.  In  computing  the  energy,  or  any  property  for  which  d>  is  an  eigen¬ 
state.  there  is  no  difference  between  these  two  averages.  This  follows  since  the 
eigenvalue  can  be  taken  out  of  the  integral  in  the  numerator  of  Fq.  -f.  In  particu¬ 
lar,  to  compute  the  energy  one  samples  the  local  energy  EL  {R  ).  Then 


<  A  >  = 


/ d>{R  )Vt(H)[Vt{!L)I{'1,t{R)}  <IR 

Jodi  )*7 .(/?  )  dR 


=  <  o  I  //  I  vk  -  > = /•;  n . 


w  hore  A  0  is  the  eigenvalue  corresponding  to  the  state  o.  The  last  equalitv 


\L 


follows  upon  noting  that  H  is  Hermitian,  and  thus  can  operate  to  the  left. 

For  expectation  values  of  quantities  whose  operators  do  not  commute  with 
H  ,  the  mixed  average  of  Eq.  5  is  only  approximate.  One  suspects  that  the  mixed 
average  is  in  some  sense  “half-way”  between  the  exact  expectation  value  (with 
respect  to  q)  and  the  variational  expectation  value,  taken  with  respect  to  the 
trial  wave  function,  i.e.  <4'r  |  .1  |  'f'y-  >.  Taken  literally,  this  implies  that 


<0  \  A  |<p>=2<4,7'  |  .4  \  <}>-<'&  t  j  .4  |  't'  7*  >  .  (7) 

This  result  can  be  readily  formalized  [19]  and  gives  an  approximate  formula  for 

expectation  values  taken  solely  with  respect  to  <p  from  just  mixed  and  variational 

averages.  However,  it  is  difficult  to  judge  how  accurate  this  approximation  is. 

Thus,  it  is  desirable  to  be  able  to  sample  exactly  from  the  distribution  jo]2. 

This  can  be  done,  though  with  some  changes  to  the  usual  Q\1C  algorithm.  The 

distribution  must  be  weighted  locally  by  6{R  )/'i>  T  (R  ).  This  quantity  is 

essentially  the  asymptotic  number  of  survivors  of  the  local  configuration  R  [20], 

Thus,  algorithmically,  one  must  follow  each  configuration  into  the  future  before 

computing  any  averages.  Our  results  have  shown  that  while  the  cariational 

approximation  is  poor,  the  approximate  formula  (7)  is  quite  good,  and  excellent 

agreement  with  exact  results  is  obtained  by  sampling  from  the  pure  |  o  j  -  1  list ri- 
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Progress  in  Current  Year 

(l)  Excited  States.  As  mentioned  in  Sect.  I.  work  with  QMC  has  been  essen¬ 
tially  limited  to  ground-state  potential-energy  surfaces  and  lowest-energy  states 
of  a  particular  symmetry  [-1-7,18].  Examples  include  our  calculation  of  the  energy 
of  the  first  excited  state  of  methylene  [7]  in  order  to  obtain  the  (until  recently) 
elusive  singlet-triplet  splitting.  Our  results  there  are  in  excellent  agreement  with 
the  most  recent  experiments.  The  restriction  on  lowest  energy  states  of  a  sym¬ 
metry  comes  from  an  essential  feature  of  the  mapping  of  the  Schrddinger  equa¬ 
tion  into  its  diffusion  analog — -namely,  that  time  in  these  two  equations  differs  by 
a  factor  of  t  .  This  means  that  the  expansion  of  a  time-dependent  molecular 
state  vector  in  energy  eigenfunctions  multiplied  by  exp (-iEt  / fi),  results  in  a  series 
in  which  only  the  lowest  energy  term  (i.e.  0)  survives  at  large  t  .  Thus  one 
obtains  exponential  convergence  to  the  lowest  energy  eigenstate. 


If  ^  t  ls  orthogonal  to  the  exact  lowest-energy  state,  one  can  see  from  Eq.  6 
that  convergence  will  be  to  the  next-lowest  energy.  (Initially  o  contains  a  super¬ 
position  of  states.)  This  fact  actually  is  used  even  in  the  calculation  of  molecular 
ground-state  energies,  as  Fermi  states  are  excited  states  (with  respect  to  the 
Boson  ground  state)  of  the  Schrodinger  equation.  By  choosing  vp  j  to  be 
antisymmetric  with  respect  to  particle  exchange,  one  can  project  out  all  sym¬ 
metric  States.  A  similar  result  holds  for  calculations  of  different  symmetry  states 
of  a  given  molecule. 
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symmetry.  This  implies  that  convergence  will  ultimately  be  to  the  lowest-energy 
state.  However,  the  fixed- node  approximation  used  to  treat  the  Fermi  problem 
may  be  used  to  advantage  here.  In  the  fixed-node  approximation,  the  nodes  of 
'I' T  are  used  to  divide  R  -space  into  separate  volume  elements.  The  Schrodinger 
equation  is  solved  separately  in  each  of  these  elements.  This  results  in  a  solution 
of  the  Schrodinger  equation  with  added  boundary  conditions.  Viewed  this  way, 
the  Fermi  problem  is  handled  by  forcing  the  generation  of  an  antisymmetric  state 
above  the  Bose  ground  state  through  the  placement  of  nodes  in  the  solution  <t>. 
In  like  manner,  other  excited  states  can  be  treated  approximately  by  imposing 
additional  nodes.  The  accuracy  of  the  approximation  will  depend  on  how  well 
these  nodes  are  placed. 

In  treating  excited  states,  traditional  ab  initio  methods  generate  wave  func¬ 
tions  with  the  correct  number  and  dimensionality  of  nodes.  Thus  such  wave 
functions  are  a  good  place  to  begin  in  searching  for  an  excited-state  ^  p  .  This 
year  we  have  used  QMC  in  this  way.  and  explored  the  type  ,,f  trial  wave  func¬ 
tions  required  for  accurate  nodal  representation  [21]. 


We  have  explored  the  effects  of  basis  set  and 
configurations  on  the  convergence  of  the  QMC  energies, 
number  of  determinants,  results  do  improve  with  the 
configurations  (see  Table  1).  Quite  possibly,  after  a  certain 
noticahle  improvements  will  result.  For  the  ground  states 
h;is  indeed  been  the  case,  where  a  single  determinant  and  a  < 
level  have  been  that  level. 


number  of  reference 
At  least  for  a  small 
number  of  reference 
level  of  calculation  no 
studied  to  date  this 
louble-zeta  expansion 
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We  are  also  exploring  a  new  approach  for  calculating  excited  states.  In  this 
method,  all  (i.e.  a  significant  number,  say  10)  excited  states  are  found  in  a  single 
calculation.  At  the  present  time,  in  order  to  avoid  problems  with  Fermi  statis¬ 
tics,  we  are  restricting  ourselves  to  vibrational  and  rotational  states  through  the 
use  of  potential  energy  surfaces.  This  restriction  can  be  lifted,  and  excited  elec¬ 
tronic  states  will  be  explored  at  a  later  time. 

The  method  involves  the  calculation  of  correlation  functions.  In  particular, 
we  calculate  the  following  matrices  of  correlation  functions 

(8) 

and 


£a/j(0=<*a(0l  l^)>  (0) 

Here  {^Q(0)}  is  a  set  of  M  trial  wave  functions,  representing  guesses  at  the  first 

M  states  of  the  system,  and  *a(t  )=c  ~r//  *a(0).  The  eigenvalues  of  E  /.V  then 

give  the  first  M  eigenvalues  of  H .  Briefly,  we  seek  the  values  of  X  that  are  solu¬ 
tions  to 

det(£-X;V)=0.  (]()) 

Fx  ponding  in  a  complete  set  of  states  {o,  }  we  have 


*  =Y 


r,'J  e  6,  . 


1  1 


I  sing  the  Binet-Cauchy  theorem  for  rectangular  matrix  products,  we  obtain 


'  1 


l: 


det(£-X,V)=  2  I  D. 

«  !<«*<.  ■•■<  I* 

where  c,  =(E, -\)e  E' '  and  |  ,u  |  =det  ]  r/'  \  .  Thus,  asymptotically  (i.e. 

/  — *oo),  the  dominant  term  in  (12)  is 


(E0-\)(El-\)...(EM-\)e-c'«  e-*"  ■  ■  ■  e —  ‘  |  D ,  M  |  ,  (13) 

and  consequently,  the  eigenvalues  are  indeed  the  first  M  eigenvalues  of  the  Ham¬ 
iltonian.  One  constraint  is  that  each  of  the  actual  M  lowest  eigenfunctions  of  H 
must  have  a  non-zero  overlap  with  at  least  one  of  the  If  this  is  not  so,  the  M 
states  obtained  will  exclude  those  states  with  zero  overlap. 

In  order  to  compute  the  correlation  matrices  of  Eqs  (8)  and  (9),  let  us  write 

E„s(<  )=<*„(£')  I  He'"  \*M)> 

=/Wr,tr)‘  vaTjW  KH  iKiR~ 

t 

-Jdl'Ei(t') 

=  <EL°{R')e  °  w\R)>  # 


-f#Enr  ) 

=  <<ELa{R’ {t +r))e  '  tif3(R{r))>  >, 


(14a) 


where  the  bracketed  quantity  [  .  .  .  ]  is  the  usual  QMC  Green’s  function  which 
propagates  the  guiding  function  (*,  )  times  the  wave  function  at  R  to 

)^(K  );  —  is  the  usual  local  energy  associated  with  a 

function  ^g  ;  wi=^/3J^g  ;  and  the  outer  average  <  .  .  .  >  in  the  final  line  is  an 
average  over  r,  corresponding  to  an  average  over  choices  of  origin  for  the  random 
walk.  For  this  last  equality  to  be  valid,  one  must  have  ergodidty  for  the  random 
walk  being  described  by  the  Green’s  function.  Replacing  Ef  in  the  final  line  by 
Kives  the  correlation  matrix  .Y,Jt  ).  Finally,  note  that  if  =  vj/  r 

this  approach  reduces  to  the  usual  QMC'. 

An  alternative  expression  that  is  perhaps  more  transparent  because  of  the 
similarity  of  its  development  to  the  usual  OMC  arguments  begins  with  Hq.  (9): 


T 


H  dR,  ,0)  ) 

E^t  )=<*„(<  )  I  H  I  )>=/■  (a  -°)  (g  ,,  )  *,  (£  « 


*,(fi  •< ) 


,  H  */0)  fjl )  __ 
<  *  (0)  * At ) 


(Mb) 


Though  different  in  form,  this  expression  is  equivalent  to  Eq.  (Ida). 


As  a  test  of  this  method,  we  have  studied  the  H90  molecule,  for  which  there 
are  a  number  of  good  potential-energy  surfaces,  as  well  as  other  theoretical  and 
experimental  results  on  vibrational  excited  states.  Our  numerical  results  are  only 
preliminary.  Using  only  2  minutes  of  VAX  8600  cpu  time,  we  found  the  ground- 
state  energy  exact  to  0.2%,  and  the  first  five  excited  states  to  within  10%.  Since 
these  are  already  the  vibrational  frequencies,  10%  is  not  too  bad.  Nevertheless, 
we  expect  these  results  can  be  significantly  improved.  The  reason  for  the  much 
larger  noise  in  the  excited  states  is  due  to  the  fact  that  all  the  trial  functions 
**(* )  approach  the  ground  state  at  large  time.  Thus,  noise  in  the  matrix  ele¬ 
ments  of  Eqs.  8  and  9  grows  exponentially. 


(2)  Molecular  Properties.  Our  work  on  properties  unrelated  to  the  energy- 
for  which  averages  from  the  true  distribution  |  £  |  2  must  be  used— has  been 
quite  successful.  The  theory  is  described  in  Sect.  1.3  above  and  will  not  be 
repeated  here.  Much  of  our  current  work  has  involved  gaining  additional  experi¬ 
ence  analogous  to  what  we  have  already  gained  in  studies  of  ground-state  ener¬ 
gies.  We  have  discovered  some  4*  T  dependence  for  certain  properties  (e.g.  dipole 
moments).  Basis-set  expansions  and  multi-determinant  y-  s  need  to  In- 
explored.  Further,  an  understanding  must  be  gained  of  why  the  averages  with 
respect  to  |  0  \  1  should  be  trial-function  sensitive.  It  is  clear  how  polarization 


wv 
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functions  can  dramatically  affect  dipole  moments  computed  with  respect  to 
!  I  2!  such  effects  should,  however,  be  irrelevant  with  the  correct  |  <p  |  2 
averaging  procedure.  The  only  approximation  is  the  nodal  positioning.  Thus  we 
also  wish  to  understand  how  placement  of  these  nodes  affects  properties.  Some 
results  for  dipole  and  quadrupole  moments  of  various  molecules  are  presented  in 
Table  2. 


(3)  Matrix  Elements  Between  Different  Electronic  States.  Thus  far  we  have 
been  discussing  the  calculation  of  expectation  values  of  the  form 
<*r  |  A  |*r>  or  <4>\  A  \'i>T>  or  <<j>  \  A  |  <£>.  Here  *r  and  0  are 
approximate  forms  for  the  same  state.  Matrix  elements  for  operators  causing 
transitions  between  different  states,  such  as  <0X  \  x  \  02> ,  are  of  particular 
importance  for  the  study  of  absorption  of  radiation  in  molecules.  We  are 
developing  QMC  methods  for  these  problems.  Theoretical  development  is 
currently  proceeding,  and  is  summarized  here.  This  effort  also  ties  in  with  our 
other  work  including  calculations  of  expectation  values  of  operators  such  as  x  . 
and  of  excited-state  properties. 


As  mentioned  above,  the  object  is  to  compute  matrix  elements  of  the  type 


do) 


rr  =  <o{  |  f(ff )  |  <p,> 

Here  /  (R  )  is  an  arbitrary  function  of  the  electronic  coordinates.  The  states  o, 
and  o.,  are  different  eigenfunctions  of  the  Hamiltonian  in  the  I WOppenheimor 
approximation.  As  described  i„  the  theory  section,  standard  QMC  applicati 
provides  expectation  values  for  a  <|uan’ti>  Q  (  /,’  )  ,,f  i|1(.  j'(,rill 


ton 


(16) 


<Q  >=fQ(R)V  Tl(Z  )4>X(R  )dR  //*  r,(£  )0i(£  )dR 

where  V j ^(K)  is  a  trial  or  guiding  function  corresponding  to  the  eigenstate  qx. 
The  idea  is  to  choose  Q  ( R  )  in  such  a  way  as  to  obtain  T . 

First  let  us  write  Q  (R_  )=h  (R_  )w  (R  ).  By  properly  choosing  the  weight 
w(R)  one  can  introduce  the  state  02.  Namely,  if  one  uses  a  trial  function  as 
a  guiding  function  to  propagate  random  walks  forward  from  R_  .  the  asymptotic 
number  of  walkers  converges  to 

w{B.)=exp[-(E2-ETi)t  ]<^rj  |  4>2>4>2{&)I*t2{R)  ■  (17) 

Here  E T  is  a  trial  energy  chosen  as  close  as  possible  to  the  energy  E 2  of  the  the 

state  <t> 2,  and  t  is  the  time  duration  that  the  walk  has  been  propagating  from 

point  &.  In  practice,  convergence  of  w  (jt_ )  is  very  rapid.  Choosing  ET~E2 

thus  leads  to 


<Q  >  = 


<*rj  ..  ..  .  *r,(£)  * 


Mr  M&m 


<*T>  I  d>!>  '*r2(£) 

It  is  thus  apparent  that  choosing  h  (R  )=f{R  )*  T  (B.  (R  )  yields 


(18) 


<Q  >  = 


<'b I  d»o> 


T  . 


(19) 


<^r,  I  <P i> 

Given  fairly  accurate,  normalized  trial  functk  ns,  the  two  overlap  integrals  in  Eq. 
(19)  should  both  be  close  to  (but  less  than)  unity.  Thus  one  can  compute  the 
desired  matrix  element.  If  necessary.  QMC  may  be  employed  to  compute  the 
overlap  integrals,  to  eliminate  this  approximation. 


Thus  far  this  method  has  only  been  tested  on  the  hydrogen  atom.  The 
results  show  that  the  method  works,  though  more  realistic  problems  need  to  be 
addressed.  Table  :j  shows  the  results  of  our  preliminary  work. 
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Table  1.  Comparison  of  the  energy  of  excited  singlet  states  of  Ho  ( B  JHU+)  and 
(£‘1E?+)  at  their  equilibrium  geometries  with  self-consistent  field  (SC'F), 
configuration  interaction  (Cl)  and  exact  results.  For  QMC,  several  trial  functions 
are  used.  They  are  constructed  from  double-zeta  (DZ)  and  double-zeta-plus- 
polarization  (DZP)  SCF  functions.  The  colun^n  headed  %  CE  gives  the  percen¬ 
tage  of  correlation  energy  recovered.  Statistical  uncertainties  are  indicated  in 
parentheses. 


Method 


MCSCFa 


U>  (U  'S*) 

*<*>  ‘'iCE  '(kcai/^iole)’ 


E(h) 


n 


il..  (E‘E/1 

C-  cp  Inexact 

(kcal/mole) 


QMC 

(DZ) 


-0.7553 


-0.748(2) 


(DZP)  -0.7536(3)  79 


8  de^*  (DZP)  -°-7563(9)  97 


Exact c 


a  Ref.  21. 
6  Ref.  22. 
f  Ref.  23. 


-0.7567 


-0.7174(8)  95 


-0.7  i  81 


Table  2.  Dipole  and  quadrupole  moments  for^  various  molecules.  Units  are  £ 

Debyes  for  the  dipole  moments  and  esu-cm2X  1CT26  for  the  quadrupole  moments.  V 

“Approximate  Formula”  refers  to  Eq.  7  of  the  text.  Statistical  uncertainties  are  i 


indicated  in  parentheses. 


Method 

Ho 

No 

Q 

Q 

C'J'T’I  A  |  ^t> 

0.49(1) 

C'J'rl  A  \4» 

0.56(2) 

Approximate 

Formula 

0.63(5) 

-1.41(20) 

<0  |  A  |  <p> 

0.61(3) 

— 

Hartree-Fock 

0.664“ 

-l.29c 

23 


Table  3.  Matrix  elements  between  the  Is  and  2s  states  of  II. 


/  (£) 

<&l3  \  f  1  <?2a  >  QMC 

<-$la  /  I  3  exact 

r 

-0.557(11) 

-0.559 

*> 

r  “ 

-2.92(7) 

-2.98 

